#include <stdio.h>
#include <unistd.h>
#include <linux/unistd.h>
#include <sys/syscall.h>
#include <sys/time.h>
#include <sys/types.h>
#include <sched.h>
#include <signal.h>
#include <setjmp.h>
#include <errno.h>
#include <assert.h>
#include <math.h>

#include "gt_include.h"

#define NUM_CPUS 4
#define NUM_GROUPS NUM_CPUS

#define NUM_THREADS 32

typedef struct matrix
{
	int m32[32][32];
	int m64[64][64];
	int m128[128][128];
	int m256[256][256];
	int rows;
	int cols;
	unsigned int reserved[2];
} matrix_t;

typedef struct __uthread_arg
{
	matrix_t *_A, *_B, *_C;
	unsigned int reserved0;

	unsigned int tid;
	unsigned int gid;
	int start_row; /* start_row -> (start_row + PER_THREAD_ROWS) */
	int start_col; /* start_col -> (start_col + PER_GROUP_COLS) */
	
}uthread_arg_t;

struct timeval tv1;
struct timeval timecreated[NUM_THREADS],timecompleted[NUM_THREADS];

static void generate_matrix(matrix_t *mat, int rows, int cols, int val)
{

	int i,j;
	mat->rows = rows;
	mat->cols = cols;
	//int **m = (int **)malloc(sizeof(int *)*mat->rows);
	//for(i = 0; i < mat->rows;i++)
	//	mat->m[i] = (int *)malloc(sizeof(int)*mat->cols);

	for(i = 0; i < mat->rows;i++)
	{
		for( j = 0; j < mat->cols; j++ )
		{
			if(rows==32)
				mat->m32[i][j] = val;
			if(rows==64)
				mat->m64[i][j] = val;
			if(rows==128)
				mat->m128[i][j] = val;
			if(rows==256)
				mat->m256[i][j] = val;
		}
	}
	return;
}

static void print_matrix(matrix_t *mat)
{
	int i, j;

	for(i=0;i<mat->rows;i++)
	{
		for(j=0;j<mat->cols;j++)
		{
			if(mat->rows==32)
				printf(" %d ",mat->m32[i][j]);
			if(mat->rows==64)
				printf(" %d ",mat->m64[i][j]);
			if(mat->rows==128)
				printf(" %d ",mat->m128[i][j]);
			if(mat->rows==256)
				printf(" %d ",mat->m256[i][j]);
		}
		printf("\n");
	}

	return;
}

static void * uthread_mulmat(void *p)
{
	int i, j, k;
	int start_row, end_row;
	int start_col, end_col;
	unsigned int cpuid;
	struct timeval tv2;

#define ptr ((uthread_arg_t *)p)

	i=0; j= 0; k=0;

	start_row = ptr->start_row;
	//end_row = (ptr->start_row + PER_THREAD_ROWS);
	end_row = (ptr->start_row + ptr->_A->rows);

//#ifdef GT_GROUP_SPLIT
//	start_col = ptr->start_col;
//	end_col = (ptr->start_col + PER_THREAD_ROWS);
//else
	start_col = 0;
	//end_col = SIZE;
	end_col = (ptr->start_col + ptr->_B->cols);
//#endif

//#ifdef GT_THREADS
	cpuid = kthread_cpu_map[kthread_apic_id()]->cpuid;
	//fprintf(stderr, "\nThread(id:%d, group:%d, cpu:%d) started",ptr->tid, ptr->gid, cpuid);
//#else
//	fprintf(stderr, "\nThread(id:%d, group:%d) started",ptr->tid, ptr->gid);
//#endif
	gt_yield();
	if(ptr->_C->rows==32)
	{
		for(i = start_row; i < end_row; i++)
			for(j = start_col; j < end_col; j++)
				for(k = 0; k < ptr->_B->rows; k++)
					ptr->_C->m32[i][j] += ptr->_A->m32[i][k] * ptr->_B->m32[k][j];
	}
	if(ptr->_C->rows==64)
	{
		for(i = start_row; i < end_row; i++)
			for(j = start_col; j < end_col; j++)
				for(k = 0; k < ptr->_B->rows; k++)
					ptr->_C->m64[i][j] += ptr->_A->m64[i][k] * ptr->_B->m64[k][j];
	}
	if(ptr->_C->rows==128)
	{
		for(i = start_row; i < end_row; i++)
			for(j = start_col; j < end_col; j++)
				for(k = 0; k < ptr->_B->rows; k++)
					ptr->_C->m128[i][j] += ptr->_A->m128[i][k] * ptr->_B->m128[k][j];
	}
	if(ptr->_C->rows==256)
	{
		for(i = start_row; i < end_row; i++)
			for(j = start_col; j < end_col; j++)
				for(k = 0; k < ptr->_B->rows; k++)
					ptr->_C->m256[i][j] += ptr->_A->m256[i][k] * ptr->_B->m256[k][j];
	}

//#ifdef GT_THREADS
	//gettimeofday(&tv2,NULL);
	//fprintf(stderr, "\nThread(id:%d, group:%d, cpu:%d) finished (TIME : %f s)",
	//		ptr->tid, ptr->gid, cpuid, (tv2.tv_sec - tv1.tv_sec + (tv2.tv_usec - tv1.tv_usec)/1000000.0));
	gettimeofday(&timecompleted[ptr->tid],NULL);
//#else
//	gettimeofday(&tv2,NULL);
//	fprintf(stderr, "\nThread(id:%d, group:%d) finished (TIME : %lu s and %lu us)",
//			ptr->tid, ptr->gid, (tv2.tv_sec - tv1.tv_sec), (tv2.tv_usec - tv1.tv_usec));
//#endif

#undef ptr
	return 0;
}

uthread_arg_t uargs[NUM_THREADS];
uthread_t utids[NUM_THREADS]; //typedef unsigned int uthread_t;

matrix_t A[NUM_THREADS], B[NUM_THREADS], C[NUM_THREADS];

void deviation_calc(double *array, int chunk_size, int chunk_num, double *mean, double *dev)
{
	int i;
	int start_index = chunk_size*chunk_num;
	int end_index = chunk_size*(chunk_num + 1);
	*mean = 0;
	*dev = 0;
	for(i=start_index; i<end_index; i++)
		*mean = *mean + array[i];
	*mean = *mean/chunk_size;
	for(i=start_index; i<end_index; i++)
		*dev = *dev + (array[i] - *mean)*(array[i] - *mean);
	*dev = sqrt(*dev/chunk_size);
	return;
}

int main()
{
	
	uthread_arg_t *uarg;
	int i, j, inx;

	gtthread_app_init();

	for(i = 0; i<NUM_THREADS; i++)
	{
		if(i>=0 && i<(NUM_THREADS/4))
		{
			generate_matrix(&A[i], 32, 32, 1);
			generate_matrix(&B[i], 32, 32, 1);
			generate_matrix(&C[i], 32, 32, 0);
		}
		if(i>=(NUM_THREADS/4) && i<(NUM_THREADS/2))
		{
			generate_matrix(&A[i], 64, 64, 1);
			generate_matrix(&B[i], 64, 64, 1);
			generate_matrix(&C[i], 64, 64, 0);
		}
		if(i>=(NUM_THREADS/2) && i<(NUM_THREADS*3/4))
		{
			generate_matrix(&A[i], 128, 128, 1);
			generate_matrix(&B[i], 128, 128, 1);
			generate_matrix(&C[i], 128, 128, 0);
		}
		if(i>=(NUM_THREADS*3/4) && i<NUM_THREADS)
		{
			generate_matrix(&A[i], 256, 256, 1);
			generate_matrix(&B[i], 256, 256, 1);
			generate_matrix(&C[i], 256, 256, 0);
		}
	}

	gettimeofday(&tv1,NULL);
	for(inx=0; inx<NUM_THREADS; inx++)
	{
		uarg = &uargs[inx];
		uarg->_A = &A[inx];
		uarg->_B = &B[inx];
		uarg->_C = &C[inx];

		uarg->tid = inx;

		uarg->gid = (inx % NUM_GROUPS);

		uarg->start_row = 0;
		uarg->start_col = 0;

		gettimeofday(&timecreated[inx],NULL);
		uthread_create(&utids[inx], uthread_mulmat, uarg, uarg->gid);
	}

	gtthread_app_exit();

	double totaltime[NUM_THREADS];
	double totaltime_mean[4], totalruntime_mean[4];
	double totaltime_dev[4], totalruntime_dev[4];
	FILE * pFile;
	pFile = fopen ("time.txt","w");
	if (pFile==NULL)
		fprintf(stderr, "Cannot create time.txt\n");
	for(inx=0; inx<NUM_THREADS; inx++)
	{
		totaltime[inx] = timecompleted[inx].tv_sec - timecreated[inx].tv_sec + (timecompleted[inx].tv_usec - timecreated[inx].tv_usec)/1000000.0;
		fprintf(pFile, "%f:%f\n", totaltime[inx], totalruntime[inx]);
	}
	fclose(pFile);
	for(i=0; i<4; i++)
	{
		deviation_calc(totaltime, NUM_THREADS/4, i, &totaltime_mean[i], &totaltime_dev[i]);
		deviation_calc(totalruntime, NUM_THREADS/4, i, &totalruntime_mean[i], &totalruntime_dev[i]);
	}

	fprintf(stderr, "TotalTime_mean");
	for(i=0; i<4; i++)
		fprintf(stderr, ":%f", totaltime_mean[i]);
	fprintf(stderr, "\n");

	fprintf(stderr, "TotalTime_dev");
	for(i=0; i<4; i++)
		fprintf(stderr, ":%f", totaltime_dev[i]);
	fprintf(stderr, "\n");

	fprintf(stderr, "CPUTime_mean");
	for(i=0; i<4; i++)
		fprintf(stderr, ":%f", totalruntime_mean[i]);
	fprintf(stderr, "\n");

	fprintf(stderr, "CPUTime_dev");
	for(i=0; i<4; i++)
		fprintf(stderr, ":%f", totalruntime_dev[i]);
	fprintf(stderr, "\n");

	return(0);
}
